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IV 


1.  SUMMARY 


Regional  seismic  phases  (e.g.,  Pn,  Sn,  Pg,  Lg )  play  a  very  important  role  in  global 
monitoring  of  low  yield  underground  nuclear  tests.  Given  a  typical  paucity  of  teleseismic 
observations  for  a  small  event,  regional  phase  arrival  times  can  be  critical  for  locating  the 
source.  Numerous  empirical  observations  have  also  shown  that  regional  phases  hold  the 
keys  to  small  event  magnitude  and  yield  estimation  and  to  discrimination  between  small 
explosions  and  earthquakes  (e.g.,  Taylor  et  al.,  1989;  Kim  et  ah,  1997;  Walter  et  ah, 
1995;  Fisk  et  al.,  1996,  2005;  Taylor,  1996;  Taylor  and  Hartse,  1997;  Hartse  et  ah,  1997; 
Patton,  2001;  Xie,  2002;  Richards  and  Kim  2007;  Zhao  et  ah,  2008).  Use  of  travel  times 
and  amplitudes  of  regional  phases  for  location,  discrimination  and  yield  estimation 
procedures  requires  that  the  regional  velocity  and  attenuation  structures  be  sufficiently 
well  approximated  (e.g.,  Myers  et  ah,  2010).  As  a  result  there  have  been  extensive  efforts 
to  determine  regional  crustal  and  upper  mantle  velocity  models  and  tomographic  models 
for  attenuation  of  regional  phases.  These  efforts  commonly  proceed  independently,  with 
geometric  spreading  models  being  used  for  simple  or  generic  structures  even  while 
detailed  regional  velocity  models  are  determined  for  travel  times.  This  situation  is 
evolving,  with  an  increasing  emphasis  on  model-based  procedures  rather  than  on 
empirical  data  trends  entering  into  operational  systems,  and  this  presents  both  challenges 
and  opportunities  for  developing  self-consistent  model-based  procedures.  The  guiding 
theme  of  this  project  is  working  toward  unified,  self-consistent  correction  of  regional 
phases  for  amplitude  and  travel  time  effects  of  heterogeneous  Earth  structure.  Our  focus 
is  on  two  of  the  critical  regional  phases,  Pn  and  Sn,  which  traverse  the  uppermost  mantle, 
and  commonly  are  assumed  to  have  very  simple  wave  behavior  controlling  their 
amplitudes  and  travel  times.  In  fact,  both  phases  have  complex  sensitivity  to  Earth 
structure  that  is  not  well  represented  by  standard  processing  assumptions  (particularly 
when  determining  attenuation  models  needed  for  the  model-based  procedures), 
motivating  use  of  regionally  appropriate  velocity  structures,  just  as  are  motivated  for 
precise  event  location.  Striving  toward  self-consistency  of  regional  structures  used  for 
travel  time  and  amplitude  analysis  is  a  logical  step  as  model-based  approaches  are 
adopted  in  the  operational  environment. 

2.  INTRODUCTION 

This  project  builds  on  a  prior  collaboration  in  which  we  developed  geometric 
spreading  corrections  for  regional  phases  Pn  and  Sn  for  simple  reference  velocity  models 
and  applied  these  to  develop  frequency-dependent  Pn  attenuation  models  for  Eurasia.  We 
confirmed  that  these  seismic  phases  have  acute  sensitivity  to  structure,  with  frequency- 
dependent  geometric  spreading  behavior  for  even  the  simplest  spherical  structures.  The 
sensitivity  to  mantle-lid  velocity  gradient  is  non-linear  and  pronounced,  but  can  be 
characterized  well  by  numerical  modeling  for  specific  structures.  Current  procedures  tend 
to  ignore  the  complex  geometric  spreading  of  these  phases  and  thus  likely  project 
artificial  frequency  dependence  into  attenuation  models  for  the  mantle  lid.  Practical 
determination  of  precise  regional  mantle-lid  velocity  gradients  is  challenging,  but 
possible  if  curvature  of  Pn  and  Sn  travel  time  branches  is  constrained  in  a  local  region. 
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This  is  exactly  what  is  being  performed  in  RSTT  (regional  seismic  travel  time) 
approaches,  although  data  scatter  and  trade-off  with  variable  crustal  legs  certainly 
reduces  confidence  in  inferred  mantle-lid  radial  velocity  gradients.  We  are  exploring  the 
predictions  of  RSTT  and  similar  models  for  amplitude  effects,  considering  both  the 
implications  for  associated  frequency-dependent  attenuation  models  and  the  prospects  of 
jointly  using  amplitude  and  travel  time  information  to  better  constrain  the  regional 
velocity  and  attenuation  structures.  Eurasian  datasets  of  regional  Pn  and  Sn  phases  will  be 
utilized,  and  our  previous  work  on  determining  frequency-dependent  geometric  spreading 
as  a  function  of  ID  velocity  model  variations  will  be  expanded  to  provide  corrections  for 
laterally  varying  regional  composites.  We  will  also  apply  a  3D  modeling  approach 
(involving  acoustic  1-way  propagators)  to  evaluate  effects  of  laterally  varying  structure 
on  the  geometric  spreading  predictions,  extending  our  current  (up  to  10  Hz)  2D  finite- 
difference  modeling  efforts  to  more  realistic  media.  Waveform  attributes,  such  as 
frequency-dependent  amplitude  ratios  of  first  Pn  arrival  to  early  Pn  coda,  that  may  help 
to  constrain  the  actual  mantle-lid  velocity  gradient  will  be  evaluated  using  the  synthetics 
and  data.  The  ultimate  goal  is  to  work  toward  procedures  for  determining  regional 
attenuation  and  velocity  models  that  are  useful  for  self-consistently  correcting  both 
travel-time  and  amplitude  measurements  for  improved  event  location,  event  identification 
and  magnitude/yield  determination. 


3.  BACKGROUND  AND  TECHNICAL  APPROACH 

Pn  and  Sn  geometric  spreading 

Accurately  accounting  for  geometric  spreading  is  critical  for  the  development  of 
meaningful  regional-phase  attenuation  models.  This  is  particularly  true  for  Pn  and  Sn 
waves  because  the  nature  of  their  wave  propagation  renders  them  acutely  sensitive  to 
uppermost  mantle  velocity  structure  and  Earth’s  sphericity.  Even  simple  one-dimensional 
(ID)  velocity  models  can  produce  geometric  spreading  of  Pn  and  Sn  that  is  strongly 
dependent  on  frequency  (e.g.,  Sereno  and  Given,  1990;  Yang  et  al.,  2007;  Avants  et  al., 
2011).  If  frequency  dependence  of  the  geometric  spreading  is  neglected,  any  inferred 
attenuation  model  based  on  using  that  geometric  spreading  will  acquire  incorrect 
frequency  dependence. 

In  our  previous  modeling  efforts,  we  generated  synthetic  seismograms  for  a  ID, 
spherical  Earth  model  (Yang  et  al.,  2007;  Avants  et  al.,  2011)  to  simulate  Pn  and  Sn 
geometric  spreading  in  such  a  structure.  We  used  the  same  generic  spherical  Earth  model 
considered  by  Sereno  and  Given  (1990)  as  the  Base  Earth  Model  for  our  simulations 
(Figure  1)  and  developed  new  Pn  and  Sn  geometric-spreading  models  using  the 
synthetics  from  the  simulation. 
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Figure  1.  Base  Earth  Model  used  for  Pn  and  Sn  simulations  and  the  development  of  new 
Pn  and  Sn  geometric-spreading  models.  Quality  factor  Q  is  infinite  throughout  the  model. 
Earth  flattening  transform  yields  the  ID  models  on  the  right,  which  have  mild  positive 
velocity  gradients  even  though  the  mantle  lid  has  constant  velocity  in  the  spherical 
structure. 

We  computed  synthetic  seismograms  using  the  source  parameterizations  and 
procedures  described  in  detail  by  Yang  et  al.  (2007)  and  extracted  the  Pn  and  Sn  portions 
of  the  synthetic  seismograms  using  fixed-velocity  windows.  We  measured  Pn  and  Sn 
amplitudes  from  the  spectra  calculated  using  the  Pn  and  Sn  synthetic  seismograms.  The 
10-Hz  Pn  amplitude  decay  in  the  spherical  Base  Earth  Model  is  shown  in  Figure  2.  Also 
plotted  in  the  figure  is  the  amplitude  decay  of  a  conical  head  wave  in  a  plane  one-layer- 
over-half-space  model  (Aki  and  Richards,  2002;  Eq.  6.26)  and  the  amplitude  decay  of 
infinite-frequency  direct  wave  in  a  spherical  Earth  model  from  ray  tracing.  At  distances 
close  to  the  critical  distance,  Pn  geometric  spreading  behaves  like  that  of  a  conical  head 
wave.  As  distance  increases,  Pn  spreading  starts  to  deviate  from  that  of  the  head  wave 
and  at  about  5°,  Pn  amplitudes  begin  to  increase.  In  the  range  between  the  critical 
distance  and  about  10°,  Pn  evolves  from  a  wave  similar  to  a  conical  head  wave  to  the 
interference  head  wave,  which  is  a  superposition  of  multiple  waves  reflected  from  the 
Moho.  At  teleseismic  distances,  the  direct-wave  spreading  approaches  that  of  the  infinite- 
frequency  diving  P  wave  from  ray  tracing. 

Pn  geometric  spreading  in  a  spherical  Earth  model  is  not  only  different  from  that  of  a 
head  wave  as  is  shown  in  Figure  2,  but  also  frequency  dependent.  Figure  3  shows  the  Pn 
amplitude-variation  surface  as  a  function  of  distance  and  frequency  for  the  simple  Base 
Earth  Model.  The  strong  frequency  dependence  of  Pn  amplitudes  is  apparent.  Amplitudes 
at  higher  frequencies  are  affected  more  by  sphericity  than  are  lower-frequency 
amplitudes. 
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epicentral  distance  (°) 


Figure  2.  10-Hz  synthetic  Pn  amplitude  decay  in  the  spherical  Base  Earth  Model  with 
constant  mantle  velocities.  The  solid  line  depicts  the  theoretical  amplitude  decay  of  a 
conical  head  wave  in  a  plane  one-layer-over-half-space  Earth  model.  The  dashed  line  is 
the  amplitude  decay  of  infinite-frequency  direct  wave  in  a  spherical  homogeneous  Earth 
model  from  ray-tracing  calculations.  From  Yang  et  al.  (2007). 


Figure  3.  Synthetic  Pn  amplitudes  as  a  function  of  epicentral  distance  and  frequency. 
From  Yang  et  al.  (2007). 
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Yang  et  al.  (2007)  presented  a  parameterized  description  of  Pn  geometric  spreading 
for  this  simple  reference  model.  The  amplitude  spectrum  of  Pn  can  be  parameterized  as 


ArAf)  =K(f)MoR(0)G{r,f)exn r  S(f) 

V  Q(f)v  J 


(1) 


with  the  new  geometric-spreading  model  expressed  as 


G(r,f)  = 


Ken 


(r0  =  1  km)  (2) 


and 


(3) 


In  Equation  1,  K  is  a  frequency-dependent  scaling  factor;  Mo  is  source  moment;  R  is 
source  radiation  pattern;  Q  is  Pn  quality  factor;  v  is  Pn  velocity;  S  is  receiver  site 
response;  r  is  epicentral  distance;  6  is  azimuth  angle  and  /  is  frequency,  ro  and  f;  are 
included  in  Equations  (2)  and  (3)  in  order  for  the  new  model  to  have  the  same  dimension 
as  standard  power-law  models.  Yang  et  al.  (2007)  derived  the  coefficients  «(/  by  fitting 
Equations  (2)  and  (3)  to  the  synthetic  data.  As  a  result,  the  new  spreading  models  better 
represent  Pn  and  Sn  geometric  spreading  in  a  more  realistic  Earth.  The  main  differences 
between  the  new  geometric-spreading  model  and  the  standard  frequency-independent 
power-law  model  are  the  addition  of  the  first  term  in  the  exponent  and  the  frequency 
dependence  of  parameters  n,. 

Pn  Spreading  Sensitivity  to  Mantle  Velocity  Gradients 

The  velocity  structure  of  the  mantle  lid  is  seldom  well  constrained.  The  most  direct 
constraint  is  provided  by  the  shape  of  the  Pn  travel  time  curve  as  a  function  of  distance 
over  a  wide-enough  distance  range  for  reliable  measurement,  and  even  when  there  is 
apparent  travel  time  branch  curvature,  it  is  difficult  to  distinguish  between  a  smooth 
gradient  with  depth  versus  constant  velocity  layering  with  small  step  increases.  Receiver 
function  methods  and  surface  wave  dispersion  inversions  sometimes  indicate  gradients  in 
the  lid  structure,  but  resolution  tends  to  be  poor,  especially  for  P-velocity.  As  a  result, 
reference  velocity  structures  often  assume  constant  or  near-constant  velocity  in  the 
mantle  lid  by  default,  so  the  Base  Earth  Model  is  not  a  bad  starting  point  in  many  cases. 
However,  there  are  regions  where  the  lid  velocity  appears  to  increase  or  decrease  with 
depth.  The  effects  of  positive  /-’-velocity  gradients  in  the  mantle  lid  on  the  frequency- 
dependent  spreading  of  Pn  are  illustrated  in  Figure  4.  The  basic  shape  of  the  spreading 
behavior  is  preserved,  but  the  distance-behavior  shifts  systematically  to  shorter  distances, 
with  overall  higher  Pn  amplitudes  as  the  gradient  increases,  ft  is  possible  to  determine 
appropriate  coefficients  for  a  spreading  representation  like  Equation  (2)  for  any  specific 
gradient,  but  this  is  only  warranted  if  there  are  a  priori  constraints  on  the  specific  model 
in  the  region  for  which  the  spreading  is  to  be  utilized. 
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Figure  4.  Pn  amplitude  versus  distance  curves  similar  to  Figure  3,  but  for  the  Base  Earth 
Model  (red  curves,  for  no  physical  gradient  in  the  lid)  and  for  two  models  with  mild 
positive  gradients  in  the  mantle  lid.  The  suite  of  curves  for  each  case  corresponds  to 
frequencies  ranging  from  0. 75  Hz  (lowest  curve)  to  12  Hz  (highest  curve).  From  Avants 
et  al.  (2011). 

Pn  Spreading  Sensitivity  to  Lateral  Mantle  Velocity  Volumetric  Heterogeneity 

In  Avants  et  al.  (2011),  we  explored  the  effects  of  lateral  heterogeneity  in  the  mantle 
lid  velocity  structure  on  Pn  spreading  using  a  2D  fourth-order  finite-difference  code  (Xie 
and  Lay,  1994).  We  included  random  lateral  volumetric  velocity  fluctuations 
characterized  by  RMS  velocity  and  varying  horizontal  and  vertical  averaging  functions 
with  different  length  scales  in  the  Model.  This  allowed  us  to  consider  2D  random 
heterogeneity  models  like  those  shown  in  Figure  5. 
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Figure  5.  Visualizations  of  2D  velocity  model  with  (top)  0.5%  (rms)  P-velocity 
heterogeneity  in  the  upper  one  hundred  kilometers  of  the  mantle  lid,  and  (bottom) 
exponential  3%  (rms)  roughness  of  the  Moho  with  40  km  scale  averaging  length.  The 
background  model  is  the  EFT  version  of  the  Base  Earth  Model,  so  there  is  a  slight 
positive  velocity  gradient  across  the  crust  and  mantle.  From  Avants  et  al.  (2011). 


The  2D  finite-different  modeling  approach  cannot  achieve  the  very  high  frequencies 
of  the  ID  wave-number  integration  method,  as  the  synthetics  are  limited  to  about  1  sec 
dominant  period.  Analysis  of  effects  of  random  heterogeneity  requires  a  statistical 
sampling  of  the  effects  associated  with  different  realizations  of  the  random  velocity 
parameters.  Many  models  have  been  run  for  various  statistical  properties  of  the  mantle  lid 
structure,  with  complete  synthetic  seismograms  being  computed. 


Figure  6  shows  the  synthetic  Pn  amplitudes  for  ensemble  averages  of  5  realizations 
each  for  variable  heterogeneity  aspect  ratios  and  rms  velocity  fluctuations  for  structures 
like  that  at  the  top  in  Figure  5,  directly  compared  to  the  result  for  the  Base  Earth  Model. 
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Figure  6.  Ensemble-averaged  Pn  amplitudes  as  a  function  of  distance  for  all 
configurations  of  random  volumetric  mantle  heterogeneity  simulated.  Pn  amplitude  for 
the  BEM  (diamonds)  are  compared  to  models  with  0.5%  RMS  (squares),  1.0%  RMS 
(triangles),  and  2%  RMS  (X's)  Vp  fluctuation  for  horizontal  (Ax)  and  vertical  (Ay) 
averaging  lengths  of  (a)  Ax  =10  km,  Ay  =10  km  (isotropic);  (b)  Ax  =  20  km,  Ay  =10 
km;  (c)  Ax  =  20  km,  Ay  =  6  km;  and  (d)  Ax  =  40  km,  Ay  =  3  km.  Each  value  is  the 
average  of  5  realizations  with  different  seed  kernels. 


An  interesting  behavior  in  Figure  6  is  that  the  presence  of  velocity  fluctuations  affects 
the  basic  shape  of  the  spreading,  and  in  some  cases  it  could  be  represented  well  by  a 
power-law  type  behavior  at  this  frequency  (~  1  Hz).  This  suggests  that  the  ID  spreading 
behavior  is  rather  delicate,  with  heterogeneity  disrupting  the  specific  interference  that 
gives  rise  to  the  complex  shape.  As  the  strength  of  the  velocity  functions  increases,  there 
is  progressive  reduction  of  curvature  of  the  amplitude-distance  trend,  and  the  overall 
amplitude  behavior  becomes  increasingly  power-law  like.  We  emphasize  that  these 
simulations  are  2D,  are  very  bandlimited  so  frequency  dependences  is  not  known,  and 
have  prescribed  velocity  heterogeneity  spectra.  Whether  they  represent,  or  even  approach 
the  prevailing  situation  in  the  mantle  lid  is  an  open  question  and  constraints  can  only 
come  from  observed  data.  In  addition,  more  accurate  modeling  will  require  broadband  2D 
and  3D  simulations.  Nevertheless,  the  results  demonstrate  the  effects  of  strong  velocity 
heterogeneity  on  Pn  geometric  spreading. 
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Pn  Spreading  Sensitivity  to  Lateral  Moho  Topography 

While  volumetric  heterogeneity  is  expected  to  result  from  geological  processes  during 
crustal  formation  and  evolution,  the  same  is  true  for  Moho  irregularities  on  multiple 
scales.  Pn  and  Sn  interactions  with  the  Moho  near  the  source  and  receiver  as  well  as 
throughout  the  propagation  and  whispering  gallery  development  provide  sensitivity  to  the 
Moho  roughness.  Avants  et  al.  (2011)  parameterized  2D  exponential  (Figure  5)  and 
Gaussian  statistical  irregularities  in  the  Moho  depth  to  explore  the  effects  of  random 
structures,  along  with  simple  step-like  structure  in  the  Moho. 

The  exponential  model  is  much  richer  in  small-scale  structure  and  produces  longer 
enduring  coda  than  the  smoother  Gaussian  model,  but  the  overall  behavior  of  the  models 
is  similar  for  the  1  sec  dominant  period  (Figure  7).  For  an  80  km  horizontal  averaging 
function  the  amplitude  decay  effects  become  significant  for  5%  RMS  heterogeneity  in  the 
depth  of  the  Moho.  There  is  only  minor  sensitivity  to  the  choice  of  horizontal  averaging 
function  for  values  larger  than  40  km  for  a  given  statistical  model. 


200  1  000  2000  200  1  000  2000  200  1000  200G 

Di sta nee  ( km)  D istan ce  (km)  D ista n ee  ( km ) 

Figure  7.  Ensemble-averaged  Pn  amplitudes  plotted  as  a  function  of  distance  for  all 
configurations  of  random  Moho  topography  heterogeneity,  grouped  by  percent  RMS 
depth  fluctuation.  Pn  amplitude  for  the  BEM  (diamonds)  are  compared  to  horizontal 
averaging  length  scales  of  20  km  (squares),  40  km  (triangles),  80  km  (X's),  and  160  km 
(asterisks),  for  percent  RMS  depth  fluctuations  of  (a)  1%;  (b)  3%;  and  (c)  5%.  From 
Avants  et  al.  (2011). 

Pn  Spreading  Model  Constrained  with  Observed  Data 

Yang  (2011)  evaluated  the  performance  of  Equation  (2)  in  correcting  observed  Pn 
amplitudes  in  Asia.  Figure  8  shows  the  1-Hz  Pn  amplitudes  after  source  and  geometric¬ 
spreading  correction  using  Equation  (2).  The  figure  shows  that  Equation  (2)  provides 
good  spreading  correction  (leaving  a  linear  trend  in  log  amplitude  versus  distance)  only 
over  a  limited  distance  range  between  about  860  and  1400  km.  Outside  this  range, 
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corrected  amplitudes  exhibit  undesirable  variations  (flattening  and  roughness)  relative  to 
a  constant  decay  rate  that  could  be  expected  for  average  attenuation  of  Pn. 
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Figure  8.  1-Hz  Pn  amplitudes  after  source  and  geometric-spreading  correction  using  the 
model  developed  by  Yang  et  al.  (2007).  The  red  line  is  the  400-point  moving  average 
showing  the  decay  trend  of  the  data. 

To  better  represent  Pn  spreading,  Yang  (2011)  developed  a  method  to  construct 
observation-based  and  region-specific  Pn  spreading  models.  In  his  method,  the  average 
Pn  Q  at  different  frequencies  is  calculated  from  source  and  spreading  (using  Equation 
(2))  corrected  Pn  amplitudes.  These  Q  values  are  then  used  to  correct  raw  Pn  amplitudes 
for  attenuation.  The  resulting  amplitude  residuals  are  then  fit  with  an  observation- 
constrained  spreading  model.  The  Pn  spreading  model  that  Yang  (2011)  developed  for 
Asia  performs  much  better,  resulting  in  a  set  of  spreading-corrected  amplitudes  that  have 
a  constant  distance  decay  slope  for  a  wide  distance  range  from  150  km  to  1600  km 
(Figure  9). 


Approved  for  public  release;  distribution  is  unlimited. 
10 


6 


-4 


-6 1 - 1 - 1 - 1 - 

200  400  600  800  1000  1200  1400  1600 

epicentral  distance  (km) 

Figure  9.  1-Hz  Pn  amplitudes  after  source  and  geometric-spreading  correction  using  the 
new  model  developed  by  Yang  (2011).  The  red  line  is  the  400-point  moving  average 
showing  the  decay  trend  of  the  data.  The  yellow  line  is  the  same  as  the  red  line  in  Figure 
8. 

The  Effect  of  Geometric  Spreading  on  Attenuation  Determination 

The  tradeoff  between  geometric-spreading  and  attenuation  makes  it  difficult  to  accurately 
estimate  one  parameter  if  the  other  parameter  is  not  well  constrained.  The  problem  is 
particularly  acute  for  Pn  and  Sn  due  to  reasons  discussed  above.  As  a  result,  an  accurate 
spreading  relationship  must  be  estimated,  even  for  individual  paths,  in  order  to  develop 
reliable  2D  attenuation  models  used  in  amplitude  corrections  for  discrimination  and  yield 
estimation. 

To  demonstrate  the  advantage  and  utility  of  a  physics-based  spreading  model  over 
traditional  power-law  model,  Figure  10  compares  the  measured  1-Hz  Pn  amplitudes  with 
the  spreading  model  of  Yang  et  al.  (2007)  and  a  power-law  model.  The  figure  shows  that 
the  Yang  et  al.  (2007)  model  has  a  consistently  slower  decay  rate  than  the  data  toward 
long  distances,  which  is  desirable  if  progressive  attenuation  occurs  to  produce  the  data 
trend  over  a  large  range  of  distance  (yellow  line).  The  power-law  model,  on  the  other 
hand,  decays  faster  than  the  data  toward  long  distances,  which  would  result  in  a  negative 
Q  when  trying  to  match  the  data  behavior. 

Figure  1 1  demonstrates  this,  showing  examples  of  Q  estimates  obtained  using  either  the 
model  of  Yang  et  al.  (2007)  or  the  power-law  model  for  multiple  frequencies.  The  model 
of  Yang  et  al.  (2007)  yields  reasonable  Q  values  at  all  frequencies  whereas  using  the 
power-law  model  results  in  either  negative  or  extremely  large  Q  values.  These  examples 
show  that  even  though  the  model  of  Yang  et  al.  (2007)  is  not  optimal,  it  provides  better 
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spreading  corrections,  due  to  its  physics  basis,  than  a  power-law  model,  which 
corresponds  to  no  known  Earth  structure. 
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Figure  10.  Comparison  of  1-Hz  Pn  amplitude  decay  with  decays  predicted  by  a  power- 
law  spreading  model  and  by  the  new  spreading  model  of  Yang  et  al.  (2007).  The  yellow 
line  is  the  400-point  moving  average  of  the  Pn  data  points  to  more-clearly  delineate  the 
amplitude-decay  trend. 

The  foregoing  review  of  accomplishments  and  enhanced  understanding  of  geometric 
spreading  for  Pn  (all  attributes  hold  for  Sn  as  well)  indicates  that  reliable  models  of 
attenuation  can  be  obtained  only  if  reliable  geometric  spreading  is  computed  for  specific 
regional  structure.  Thus,  a  strategy  for  determining  regional  structures  with 
observationally-constrained  mantle-lid  velocity  gradients  is  needed. 

The  strategy  of  regional  seismic  travel  time  (RSTT)  model  parameterization  advanced 
by  Myers  et  al.  (2010)  adopts  a  geographic  node  description  of  crustal  and  upper  mantle 
models,  where  the  mantle  portion  is  parameterized  by  a  linear  P-wave  velocity  gradient. 
This  parameterization  is  very  similar  to  that  of  the  Base  Earth  Model  used  in  the  Pn  and 
Sn  geometric  spreading  calculations  of  Yang  et  al.  (2007)  and  Avants  et  al.  (2011)  (Figure 
12).  Myers  et  al.  (2010)  used  regional  Pn  travel  times  to  develop  a  tomographic  model  of 
Pn  velocity  variation  and  mantle-gradient  variation  in  the  mantle  lid  (Figure  13).  While 
the  parameterization  is  not  explicitly  intended  to  predict  amplitude  behavior,  nor  are  the 
resulting  velocity  models  believed  to  be  accurate  in  each  case,  the  models  can  be  used  as 
a  starting  point  to  explore  amplitude  behavior. 
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Figure  1 1 .  Comparison  of  average  Q  estimates  from  Pn  amplitudes  corrected  for  source 
and  geometric-spreading  effects  using  different  spreading  models.  Q  estimates  at  1,  4  and 
8  Hz  are  compared. 
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Figure  12.  Global  model  parameterization,  (a)  An  example  tessellation  with 
approximately  5°  grid  spacing.  The  inset  shows  the  1°  used  in  the  study  by  Myers  et  al. 
(2010).  Color  indicates  Moho  depth  of  the  starting  model,  (b)  An  example  velocity  vs. 
depth  profile  as  defined  at  each  node.  The  mantle  portion  of  the  profile  is  specified  by  the 
velocity  at  the  crust/mantle  interface  and  a  linear  gradient.  From  Myers  et  al.  (2010). 

Regional  phase  arrival  times  were  used  by  Myers  et  al.  (2010)  to  develop  a 
tomographic  model  (Figure  13)  for  variation  in  velocity  of  the  uppermost  mantle  and 
variation  in  the  strength  of  the  mantle  gradient.  This  provides  regional  models  that  fit 
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travel  times,  but  the  precision  of  the  velocity  gradients  in  the  mantle  is  impacted  by 
scatter  in  the  data,  averaging  of  nodes  in  a  given  region,  and  lack  of  resolution  of 
curvature  of  Pn  travel  times.  While  the  parameterization  is  not  explicitly  intended  to 
predict  amplitude  behavior,  the  models  can  be  used  to  do  so. 

We  use  RSTT  derived  models  and  associated  travel  time  data,  with  corresponding 
regional  phase  amplitude  measures  as  a  function  of  frequency  to  evaluate  use  of  travel 
time  models  as  first-order  structures  for  computing  regional  geometric  spreading 
corrections  (thereby  having  unified  model  frameworks  for  travel  time  and  amplitudes), 
and  as  a  starting  point  for  revision  of  the  models  guided  by  observed  amplitude  behavior 
that  is  not  well  matched. 
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Figure  13.  Comparison  of  starting  and  RSTT  models,  (a)  Velocity  below  the  Moho  for 
starting  model  and  (b)  RSTT  model,  (c)  Mantle  gradient  (km/sec/km)  for  starting  model 
and  (d)  RSTT  model.  From  Myers  et  al.  (2010). 
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4.  RESULTS  AND  DISCUSSION 


Extending  the  bandwidth  of  2D  finite-difference  calculations  for  large-distance  Pn 
synthetics  was  one  of  the  initial  priorities.  The  fourth-order  code  used  by  Avants  et  al. 
(2011)  was  adapted  to  and  tested  on  our  cluster,  and  stable  P-SV  and  SH  calculations 
were  achieved  for  ranges  of  up  to  1000  km  and  frequencies  up  to  10  Hz,  for  models  with 
a  flat  free-surface  and  laterally  varying  structure,  either  with  or  without  small-scale 
statistical  heterogeneity  being  included.  This  enables  examination  of  Pn  geometric 
spreading  behavior  in  laterally  varying  2D  structure  over  broader  bandwidth  than  in  prior 
work.  Figure  14  shows  wavefield  snapshots  for  simple  ID  structures  out  to  ranges  of  900 
km.  These  models  have  either  constant  velocity  or  linear  gradients  in  the  mantle  lid. 
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Figure  14.  Snapshots  of  wavefields  computed  for  ID  structures  (constant  lid  gradient  is 
w/o  EFT)  using  the  enhanced  2D  finite-difference  code  implementation. 
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Corresponding  synthetic  seismogram  profiles  are  shown  in  Figure  15,  with  the 
amplitudes  being  clipped  for  the  very  large  amplitude  direct  Pg  arrivals  at  closer 
distances  to  allow  the  Pn  signal  to  be  seen  clearly.  The  reduction  velocity  is  8.1  km/s, 
equal  to  the  constant  mantle  velocity  layer  in  the  BEM  case  with  no  earth-flattening. 


(a)  Const  lid  w/o  EFT 


(b)  Const  lid 


(c)  Gradient  0.001 


(d)  Gradient  0.002 


Figure  15.  Vertical  component  synthetic  seismograms  computed  by  the  2D  finite- 
difference  code  plotted  with  a  reduction  velocity  of  8.1  km/s,  which  is  the  constant  lid 
velocity  in  the  basic  Earth  model  (BEM).  Note  the  small  waveform  changes  associated 
with  positive  lid  velocity  gradients  and  the  systematic  travel  time  curvature  that  would 
allow  the  velocity  gradient  to  be  resolved  given  sufficient  regional  observations. 
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With  the  finite-difference  code  being  fully  tested  and  artifacts  from  grid  dispersion  being 
suppressed  by  very  fine-model  representation  of  the  medium,  we  computed  synthetics  for 
ID  reference  models  (Figure  16)  and  2D  structures  with  slowly  varying  lateral  changes  in 
lid  gradient.  The  goal  was  to  establish  the  behavior  of  the  frequency  dependent  Pn 

geometric  spreading  for  simple  models  like  those  characterized  by  RSTT. 

Vp 


Vp  (km/s) 

5.5  6.0  6.5  7.0  7.5  8.0  8.5  9.0  9.5 


Figure  16.  Reference  velocity  models  considered  as  ID  structures  or  end-members  of  2D 
structures  with  laterally  varying  mantle  lid  gradients  for  broadband  Pn  waveform 
computations  up  to  10  Hz  at  1000  km  distance.  At  depth  the  models  tie  to  the  IASP-91 
structure  used  as  a  reference  in  most  tomography  or  location  procedures. 
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The  travel  time  variations  for  these  models  are  small,  but  potentially  measureable,  with 
picks  from  the  synthetic  waveforms  being  plotted  versus  distance  in  Figure  17. 
Potentially,  regional  observations  could  suffice  to  constrain  the  lid  gradient,  if  the 
structure  is  laterally  uniform,  if  crustal  and  site  heterogeneity  is  accounted  for,  and  if 
enough  data  are  available  to  detect  the  subtle  curvature  of  the  travel  time  curves  (versus 
simply  mapping  a  different  constant  lid  velocity  value  that  would  take  out  the  linear  trend 
in  the  time  picks. 
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Figure  17.  Arrival  time  picks  from  profiles  of  synthetics  in  Figure  15  displaying  the  weak 
curvature  introduced  by  linear  mantle  lid  gradients.  While  the  variations  are  subtle,  given 
sufficient  observations  and  laterally  uniform  structure,  the  mantle  gradient  can  be 
constrained  by  Pn  arrival  time  picks,  providing  velocity  models  that  can  be  used  for  self- 
consistent  regional  amplitude  spreading  calculations. 


For  the  ID  models  shown  in  Figure  16,  we  computed  frequency-dependent  geometric 
spreading  amplitude  decay  out  to  1000  km  using  the  2D  model,  correcting  for  behavior  of 
a  point-source  versus  a  line-source.  This  provides  the  amplitude  decay  curves  shown  in 
Figure  18.  These  patterns  are  very  similar  to  wavenumber-integration  calculations  for  the 
same  ID  velocity  models  computed  in  Avants  et  al.  (2011)  shown  in  Figure  4.  The 
complex  frequency-dependence  of  Pn,  with  minima  in  the  amplitude-distance  curves  and 
rising  amplitudes  at  larger  distances,  particularly  for  higher  frequencies  being  evident  for 
all  structures  with  a  positive  gradient  in  the  ID  flat  earth  structure.  For  lid  gradients  in 
excess  of  0.002  sec'1,  the  distance  and  frequency  dependence  are  quite  dramatic  as  shown 
before,  with  complex  plateauing  and  oscillation  of  the  high  frequency  components.  Some 
of  the  large  distance,  high  frequency  behavior  is  not  well-resolved,  but  very  similar 
behavior  is  found  for  the  completely  independent  frequency-wavenumber  integration 
synthetics  in  Figure  4.  Overall,  these  calculations  provide  a  good  level  of  confidence  in 
the  broadband  Pn  synthetics,  and  we  proceed  to  consider  results  for  2D  models. 
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Figure  18.  Pn  amplitude  decay  as  a  function  of  narrowband  frequency  versus  distance  for 
models  in  Figure  16. 
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Figure  19.  Pn  amplitude  decay  as  a  function  of  narrowband  frequency  versus  distance  for 
models  in  Figure  16. 
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As  an  initial  exploration  of  frequency-dependence  of  Pn  spreading  in  laterally  varying 
structures,  we  computed  2D  synthetics  for  models  where  the  mantle  lid  velocity  gradient 
varied  laterally.  Example  structures  are  shown  in  Figure  20.  These  approximate  regions 
where  average  ID  structures  differ  beneath  source  and  receiver  regions.  We  computed 
cases  for  sources  and  receivers  interchanged  to  explore  asymmetry  in  the  Pn  behavior. 

Distance  (km) 


0  250  500  750  100C 
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Figure  20.  2D  velocity  models  with  laterally  varying  mantle  lid  velocity  gradients  that  vary 
beneath  the  source  and  receiver  locations.  Finite-difference  Pn  synthetics  were  computed  out 
to  ranges  of 1000  km  and  for  frequencies  up  to  10  Hz. 
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Normalized  amplitude  Normalized  amplitude  Normalized  amplitude  Normalized  amplitude 


Amplitude  distance  curves  are  displayed  for  varying  frequencies  in  Figure  21,  showing  the 
differences  from  ID  models  with  zero  or  constant  velocity  gradients  and  the  2D  cases  that  have 
the  source  or  receiver  in  regions  with  the  two  end-member  ID  models. 


(a)  0.000  (b)  0.000-0.001 
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(c)  0.002-0.000  (d)  0.002  —  4.0 -6.0  Hz 


Figure  21.  Pn  amplitude  versus  distance  for  ID  (Figure  16)  and  2D  models  (Figure  20) 
for  varying  bandpass  center  frequencies. 
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The  clear  tendency  is  for  the  Pn  spreading  along  the  path  to  be  dominated  by  the  source 
side  structure  at  closer  distances  and  the  receiver  side  structure  at  longer  distances,  with 
overall  pattern  being  very  different  when  the  direction  on  the  path  is  reversed.  This 
behavior  places  special  emphasis  on  working  toward  path-specific  corrections,  as  any 
regional  average  used  for  all  paths  in  the  region  would  clearly  provide  a  poor 
approximation.  The  path  averaging  behavior  is  further  emphasized  in  Figure  22,  which 
shows  the  results  for  just  a  single  frequency  (3.75  Hz).  Note  the  differences  between  the 
red  and  green  curves,  which  are  for  simple  interchange  of  the  source  and  receiver 
positions.  Further  exploration  of  laterally  varying  models  and  addition  of  small-scale 
scattering  structure  will  be  pursued  along  with  exploring  3D  calculations  using  a  one-way 
propagator  approach. 


Distance  (km) 
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0.000 


Gradient  0.000 
Gradient  0.000  -  0.002 
Gradient  0.002  -  0.000 
Gradient  0.002 


Figure  22.  Pn  amplitude  versus  distance  for  ID  (Figure  16)  and  2D  models  (Figure  20) 
for  a  3.75  Hz  center  frequency  of  a  narrow  bandpass.  The  near-source  velocity  gradient 
is  indicated  first  in  each  laterally  varying  calculation. 
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In  parallel  with  the  computational  effort,  we  have  empirically  examined  the  use  of 
regionalized  corrections  for  Pn  amplitudes  guided  by  RSTT  regional  velocity  models  and 
by  empirical  analysis  of  Pn  data  with  the  constraint  of  realistic  spreading  behavior  as 
established  by  Yang  (2011).  Figure  23  displays  the  regionalization  guided  by  RSTT 
models. 


x  10 


Figure  23.  RSTT  model  lid  velocity  gradients  (top)  and  average  Moho  P  velocity 
(bottom)  with  a  superimposed  grid  developed  to  define  provinces  of  similar  parameters. 
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The  LANL  database  for  Pn  measurements  (Figure  24)  is  used  to  extract  signals  and 
measured  Pn  amplitudes  in  varying  passbands  to  define  the  sampling  of  regions  (Figure 
25)  for  which  RSTT  mantle  lid  gradient  and/or  average  Moho  velocity  are  relatively 
uniform. 


Figure  24.  LANL  Eurasian  Pn  data  base  coverage  showing  Pn  paths  sampling  the  RSTT- 
based  regionalization. 


HJP  provinces  and  regions 


Figure  25.  Regional  numbering  scheme  for  the  RSTT-based  subdivision.  Pn 
observations  in  each  region  have  been  evaluated.  Model  spreading  parameters  for  each 
RSST  structure  have  been  determined. 
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For  each  region  identified  in  Figure  25,  we  computed  spreading  behavior  for  a  ID 
average  of  the  local  RSTT  model  using  wavenumber  integration,  and  applied 
corresponding  frequency-dependent  geometric  spreading  corrections  to  the  MDAC- 
corrected  Pn  observations  for  that  region.  We  also  applied  simple  power-law  spreading 
corrections  for  a  power  of  -1.1,  and  we  determined  region-specific  spreading  corrections 
following  the  procedure  of  Yang  (2011).  This  gives  three  sets  of  corrected  Pn  data  for 
each  region.  Plotting  these  sets  of  corrected  Pn  amplitudes  (log  10)  as  functions  of 
distance,  allows  an  estimation  of  the  average  Q  for  that  region  for  each  frequency.  We 
display  results  for  regions  4  and  13  (Figure  25)  in  Figures  26a, b  and  27a, b,  respectively. 
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Region  4:  RSTT  spreading  correction 


Region  4:  Power-law  spreading  correction 


Region  4:  RSTT  spreading  correction 


Region  4:  Region-specific  spreading  correction 
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Region  4:  Region-specific  spreading  correction 


Figure  26a.  Pn  amplitudes  versus  distance  for  Region  4  with  corrections  using  a  power- 
law  (top  row),  RSTT-based  synthetic  computations  (middle  row)  or  empirical  region- 
specific  corrections  following  the  procedure  of  Yang  (2011)  (bottom  row)  for  passband 
center  frequencies  of  1.06  Hz  (left)  and  4.24  Hz  (right).  The  measured  average  Q  for 
each  case  obtained  from  the  regression  curve  is  indicated  in  the  upper  right. 

Each  set  of  geometric  spreading  corrections  results  in  relatively  smooth  linear  amplitude 
versus  distance  trends  that  suggest  a  corresponding  average  Q  for  that  frequency.  But  the 
Q  values  that  are  inferred  vary  dramatically,  depending  on  the  geometric  spreading  used. 
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Region  4:  Power- law  spreading  correction 


Region  4:  Region-specific  spreading  correction 
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Figure  26b.  Pn  amplitudes  versus  distance  for  Region  4  (map  inset)  with  corrections 
using  a  power-law  (top  row),  RSTT-based  synthetic  computations  (middle  row)  or 
empirical  region-specific  corrections  following  the  procedure  of  Yang  (2011)  (bottom 
row)  for  passband  center  frequency  of  8.94  Hz  (right).  The  measured  average  Q  for  each 
case  obtained  from  the  regression  curve  is  indicated  in  the  upper  right. 

In  general,  as  found  by  Yang  et  al.  (2007)  and  Yang  (2011),  use  of  a  ‘standard’  power- 
law  spreading  correction  with  an  exponent  of  -1.1  that  does  not  vary  with  frequency  (and 
for  which  there  is  no  specific  known  structure)  gives  high  Q  values  in  every  case,  and 
even  negative  Q  value  for  Region  13  at  1.06  Hz.  For  Region  4,  the  RSTT  velocity 
gradient  is  relatively  high  (0.0026  s'1)  while  for  Region  13  it’s  average  (0.0012  s'1).  In 
both  cases,  the  addition  of  the  Earth  flattening  transformation  (EFT)  gradient  of  0.0013  s' 
1  gives  positive  velocity  gradients  that  produce  strong  frequency-dependent  Pn  amplitude 
corrections.  The  data  have  somewhat  less  scatter  and  systematically  lower  Q  values  are 
inferred  for  each  region  using  the  RSTT  based  corrections.  Q  values  that  are  1.6- 1.8 
times  higher  than  the  RSTT-based  results  are  found  for  corrections  based  on  the  region- 
specific  spreading  corrections  following  Yang  (2011).  The  absolute  values  and  the 
frequency-dependence  of  the  attenuation  parameters  are  more  similar  for  the  RSTT-  and 
regional-based  spreading  corrections.  For  these  regions  there  is  no  clear  reduction  in 
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scatter  to  favor  a  specific  spreading  model,  but  the  Q  values  appear  to  be  more  reasonable 
for  either  of  the  frequency-dependent  spreading  corrections  results.  As  long  as  the  Q 
model  for  a  given  spreading  model  is  used,  total  amplitude  calculations  will  be 
reasonably  similar.  The  main  advantages  of  using  a  more  appropriate  local  velocity 
model  to  process  the  amplitudes  is  that  a  more  reliable  Q  model  may  be  estimated  and 
that  could  be  valuable  for  application  to  full  waveform  modeling,  for  extrapolation  to 
poorly  sampled  paths,  and  for  overall  self-consistency  of  event  locations  and  amplitudes. 


Region  13:  Power-law  spreading  correction 
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Region  13:  RSTT  spreading  correction 
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Region  13:  RSTT  spreading  correction 


Region  13:  Region-specific  spreading  correction 


Figure  27a.  Pn  amplitudes  versus  distance  for  Region  13  with  corrections  using  a  power- 
law  (top  row),  RSTT-based  synthetic  computations  (middle  row)  or  empirical  region- 
specific  corrections  following  the  procedure  of  Yang  (2011)  (bottom  row)  for  passband 
center  frequencies  of  1 .06  Hz  (left)  and  4.24  Hz  (right).  The  measured  average  Q  for 
each  case  obtained  from  the  regression  curve  is  indicated  in  the  upper  right. 
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Region  13:  Power-law  spreading  correction 


Region  13:  Region-specific  spreading  correction 


Figure  27b.  Pn  amplitudes  versus  distance  for  Region  13  (map  inset)  with  corrections 
using  a  power-law  (top  row),  RSTT-based  synthetic  computations  (middle  row)  or 
empirical  region-specific  corrections  following  the  procedure  of  Yang  (2011)  (bottom 
row)  for  passband  center  frequency  of  8.94  Hz  (right).  The  measured  average  Q  for  each 
case  obtained  from  the  regression  curve  is  indicated  in  the  upper  right. 
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5.  PRELIMINARY  FINDINGS 

New  2D  computations  of  Pn  geometric  spreading  have  established  that  the  amplitude 
variation  is  dependent  on  direction  of  the  path  and  that  complex  frequency  dependence  is 
expected  to  have  mixing  of  the  near-source  and  near-receiver  structures.  This  further 
weakens  the  rationale  for  using  average  power-law  representations  of  the  Pn  (and  Sn) 
geometric  spreading.  As  the  only  structure  that  will  result  in  a  power  law  behavior  for 
sure  are  those  with  negative  critical  gradients  in  the  mantle  lid  (a  possible,  but  likely  very 
rare  situation),  we  further  support  Yang’s  (2011)  contention  that  a  simple  velocity  model 
based  correction,  or  an  empirically  grounded  correction  is  a  more  sensible  choice  than  a 
frequency- independent  power  law.  We  will  explore  some  models  with  fine-scale 
heterogeneity  to  evaluate  whether  the  tendency  to  smooth  out  the  distance  behavior  for  1 
Hz  signals  as  found  by  Avants  et  al.  (2011)  holds  for  broadband  Pn  signals. 

Application  of  regionalized  Pn  spreading  corrections  based  on  RSTT  models  (derived 
exclusively  from  regional  travel  time  observations  and  velocity  models)  or  regionalized 
fitting  of  attenuation  models  to  data  provide  fairly  similar  Q(f)  estimates,  much  lower 
than  those  from  constant  power  law  spreading  coefficient.  While  the  choice  of  optimal 
spreading  is  not  unambiguous  in  every  case  based  on  some  measure  of  variance  reduction 
in  the  model  performance,  the  more  reasonable  attenuation  parameters  obtained  using  the 
RSTT  or  empirical  spreading  corrections  suggest  they  are  more  reasonable  than  standard 
power  law  assumptions. 
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